       subroutine ydrov (li,lf,dt,wg,w2,wr,pl,veg,wrmax,xl,tcdil
     :   ,lev,letr,leg,ruigc,ruitc,vsw,cw1,cw2,wgeq,xdd2)
c
      implicit double precision (a-h,o-z)
c
      real*8 leg,lev,letr
      dimension wg(2),w2(2),wr(2)
c
c       evolution of wr
c
!      write(*,*) 'ydrov>>',li,lf,dt,wg(li),w2(li),wr(li)
!    :   ,pl,veg,wrmax,xl,tcdil
!    :   ,lev,letr,leg,ruigc,ruitc,vsw,cw1,cw2,wgeq,xdd2

       zev = lev / xl
       zetr = letr / xl


       zeg = leg / xl
       wr(lf) = wr(li) - dt * (zev - zetr - veg * pl )
       wr(lf) =dmax1(0.d0, wr(lf))
       ztemp= (wr(lf) - wrmax) / dt
       ruir =dmax1(0.d0, ztemp)
       wr(lf) = min(wr(lf), wrmax)
       zpg = (1. - veg) * pl + ruir
c
c         evolution of wg
c       - - - - - - - - - - - - - -
c
c      write(*,*) 'ydrov:',wg(li),tcdil,cw2,wgeq,dt,cw1,zeg,zpg
       wg(lf) = (wg(li) -
     &    dt * (cw1 * (zeg - zpg) / 1000. - cw2 * wgeq / tcdil))
     &    /(1.+dt*cw2/tcdil)
c
       ztemp=wg(lf)-vsw
       ruig =dmax1(0.d0, ztemp)
       ruigc = ruigc + ruig * 10.

c        evolutlon of w2
c       _______________
c
      ruip = ruig * 10./(xdd2*1000.)
      w2(lf) = w2(li) - dt*(zeg + zetr - zpg + ruip) / (xdd2 * 1000.)

c     ---------------------------------------------------
c                Run-off calculations
c     ----------------------------------------------------

      ztemp = w2(lf) - vsw
      ruitc = ruitc +dmax1(0.d0,ztemp) * xdd2 * 1000. + ruig * 10.
      wg(lf) = min(wg(lf), vsw)
      w2(lf) = min(w2(lf), vsw)

      return
      end